home *** CD-ROM | disk | FTP | other *** search
/ Graphics Plus / Graphics Plus.iso / general / modelers / geomview / source.lha / Geomview / src / lib / gprim / mesh / meshnormal.c < prev    next >
Encoding:
C/C++ Source or Header  |  1993-02-08  |  4.6 KB  |  177 lines
  1. /* Copyright (c) 1992 The Geometry Center; University of Minnesota
  2.    1300 South Second Street;  Minneapolis, MN  55454, USA;
  3.    
  4. This file is part of geomview/OOGL. geomview/OOGL is free software;
  5. you can redistribute it and/or modify it only under the terms given in
  6. the file COPYING, which you should have received along with this file.
  7. This and other related software may be obtained via anonymous ftp from
  8. geom.umn.edu; email: software@geom.umn.edu. */
  9.  
  10. /* Authors: Charlie Gunn, Stuart Levy, Tamara Munzner, Mark Phillips */
  11.  
  12. #include "meshP.h"
  13. static void mnorm();
  14.  
  15. Mesh *
  16. MeshComputeNormals(Mesh *m)
  17. {
  18.     register HPoint3 *normp = NULL, *normptr, *pptr;
  19.     int i;
  20.     if(m->n) {
  21.         GeomFree(m->n);
  22.     }
  23.     m->n = GeomNewN(Point3, m->nu * m->nv);
  24.     m->flag |= MESH_N;
  25.     /* if this a '4D' mesh, we need to dehomogenize (temporarily)
  26.     before computing the normals */
  27.     if (m->geomflags & VERT_4D)    {
  28.         normp = GeomNewN(HPoint3, m->nu * m->nv);
  29.         for (normptr = normp, pptr = m->p, i=0; 
  30.         i < m->nu * m->nv; 
  31.         ++i, ++pptr, ++normptr) 
  32.             HPt3Normalize(pptr, normptr);
  33.         }
  34.     mnorm(m->geomflags & VERT_4D ? normp : m->p, m->n, m->nu, m->nv,
  35.          m->flag & MESH_UWRAP, m->flag & MESH_VWRAP, m->flag & MESH_EVERT);
  36.     if (m->geomflags & VERT_4D) GeomFree(normp);
  37.     return m;
  38. }
  39.  
  40.  
  41. static void
  42. mnorm(ap, an, nu, nv, uwrap, vwrap, evert)
  43.     HPoint3    *ap;
  44.     Point3    *an;
  45.     int    nu, nv;
  46.     int    uwrap, vwrap;
  47.     int    evert;
  48. {
  49.     register HPoint3 *prev, *next;
  50.     register Point3 *n;
  51.     register int k;
  52.     int u, v;
  53.     float x,y,z, norm;
  54.     float unit;
  55.  
  56.     /*
  57.      * We set the normal at each point to be the mean of the
  58.      * cross products at the four adjacent points.  I.e.:
  59.      *   n        Given mesh point p and its four orthogonal neighbors
  60.      * w p e    e, n, w, s we define displacements E=e-p, N=n-p, ...
  61.      *   s        and compute normal(p) = ExN + NxW + WxS + SxE
  62.      * This turns out to be equal to (e-w) x (n-s) -- independent of p!
  63.      *
  64.      * Since the index arithmetic gets a bit messy at the boundaries
  65.      * we make two passes.  Pass 1 computes e-w and holds it in
  66.      * the array to be filled with normals.
  67.      * Pass 2 computes n-s, takes the cross product and normalizes.
  68.      */
  69.  
  70.     unit = evert ? -1.0 : 1.0;
  71.  
  72.     for(u = 0; u < nu; u++) {
  73.         if(u == 0) {
  74.         prev = &ap[uwrap ? nu-1 : 0];
  75.         next = &ap[u+1];
  76.         } else if(u == nu-1) {
  77.         prev = &ap[u-1];
  78.         next = &ap[uwrap ? 0 : u];
  79.         } else {
  80.         prev = &ap[u-1];
  81.         next = &ap[u+1];
  82.         }
  83.         n = &an[u];
  84.         k = nv;
  85.         do {
  86.         n->x = next->x - prev->x;    /* e - w */
  87.         n->y = next->y - prev->y;
  88.         n->z = next->z - prev->z;
  89.         n += nu;            /* advance to next v row */
  90.         prev += nu;
  91.         next += nu;
  92.         } while(--k > 0);
  93.     }
  94.  
  95.     for(v = 0; v < nv; v++) {
  96.         if(v == 0) {
  97.         prev = &ap[vwrap ? nu*(nv-1) : 0];
  98.         next = &ap[nu*1];
  99.         } else if(v == nv-1) {
  100.         prev = &ap[nu*(v-1)];
  101.         next = &ap[vwrap ? 0 : nu*v];
  102.         } else {
  103.         prev = &ap[nu*(v-1)];
  104.         next = &ap[nu*(v+1)];
  105.         }
  106.         n = &an[nu*v];
  107.         k = nu;
  108.         do {
  109.         Point3 t;
  110.  
  111.         t.x = next->x - prev->x;
  112.         t.y = next->y - prev->y;
  113.         t.z = next->z - prev->z;
  114.         x = n->y * t.z - n->z * t.y;
  115.         y = n->z * t.x - n->x * t.z;
  116.         z = n->x * t.y - n->y * t.x;
  117.         norm = x*x + y*y + z*z;
  118.         if(norm == 0.0) {
  119.             /*
  120.              * Oh no, degenerate norm.
  121.              * Let's hope it happened because (part of) a row of
  122.              * mesh points coincided -- maybe even though N=S,
  123.              * NE != SE or NW != SW.
  124.              */
  125.             if(t.x == 0.0 && t.y == 0.0 && t.z == 0.0) {
  126.             if(k > 1) {
  127.                t.x = (next+1)->x - (prev+1)->x;
  128.                t.y = (next+1)->y - (prev+1)->y;
  129.                t.z = (next+1)->z - (prev+1)->z;
  130.             }
  131.             if(t.x == 0.0 && t.y == 0.0 && t.z == 0.0 && k < nu) {
  132.                t.x = (next-1)->x - (prev-1)->x;
  133.                t.y = (next-1)->y - (prev-1)->y;
  134.                t.z = (next-1)->z - (prev-1)->z;
  135.             }
  136.             }
  137.             if(n->x == 0.0 && n->y == 0.0 && n->z == 0.0) {
  138.             /* Do likewise in E-W direction. */
  139.             register HPoint3 *cur = &ap[nu*(v+1) - k];
  140.  
  141.             if(k == 1) cur--;
  142.             else if(k == nu) cur++;
  143.             if(v > 0) {
  144.                 cur -= nu;        /* SE-SW */
  145.                 n->x = (cur+1)->x - (cur-1)->x;
  146.                 n->y = (cur+1)->y - (cur-1)->y;
  147.                 n->z = (cur+1)->z - (cur-1)->z;
  148.             }
  149.             if(n->x==0.0 && n->y==0.0 && n->z==0.0 && v < nv-1) {
  150.                 cur += 2*nu;    /* NE-NW */
  151.                 n->x = (cur+1)->x - (cur-1)->x;
  152.                 n->y = (cur+1)->y - (cur-1)->y;
  153.                 n->z = (cur+1)->z - (cur-1)->z;
  154.             }
  155.             }
  156.             x = n->y * t.z - n->z * t.y;
  157.             y = n->z * t.x - n->x * t.z;
  158.             z = n->x * t.y - n->y * t.x;
  159.             norm = x*x + y*y + z*z;
  160.             if(norm == 0.0) {
  161.             /* Oh well. */
  162.             n->x = unit;
  163.             norm = 1.0;
  164.             }
  165.         }
  166.         norm = unit / sqrt(norm);
  167.         n->x = x*norm;
  168.         n->y = y*norm;
  169.         n->z = z*norm;
  170.  
  171.         n++;            /* Next u column */
  172.         prev++;
  173.         next++;
  174.         } while(--k > 0);
  175.     }
  176. }
  177.